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Abstract 

We study the evolution of gravitational waves(GWs) after inflation in a brane-world 
cosmology embedded in five-dimensional anti-de Sitter spacetime. Contrary to the 
standard four-dimensional results, the GWs at the high-energy regime in brane- 
world model suffer from the effects of the non-standard cosmological expansion 
and the excitation of the Kaluza-Klein modes(KK-modes), which can affect the 
amplitude of stochastic gravitational wave background significantly. To investigate 
these two high-energy effects quantitatively, we numerically solve the wave equation 
of the GWs in the radiation dominated epoch at relatively low-energy scales. We 
show that the resultant GWs are suppressed by the excitation of the KK modes. 
The created KK modes are rather soft and escape away from the brane to the bulk 
gravitational field. The results are also compared to the semi-analytic prediction 
from the low-energy approximation and the evolved amplitude of GWs on the brane 
reasonably matches the numerical simulations. 
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1 Introduction 



The stochastic gravitational wave background (stochastic GWB) generated 
during the inflationary era is a promising source for the low-frequency gravi- 
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tational waves and this can provide a direct way to probe the extremely early 
Universe. In particular, such GWs are expected to have information about the 
extra dimensions. Inspired by the recent development in particle physics, there 
has been a lot of active debate on the possibility that we live in a brane which is 
a three dimensional world embedded in a higher dimensional space [1][2][3][4]. 
If this is true, the standard prediction of stochastic GWB from the four- 
dimensional theory is dramatically altered and the GWB can be a powerful 
tool to discriminate the presence/absence of the extra-dimension. For instance, 
in the case of the single-brane model proposed by Randall & Sundrum[2], the 
scale of the extra dimension is well below the length / < 0.1mm according to 
the current experiments of the Newton gravity (e.g., [5]). Therefore, the effect 
of the extra-dimension can be imprinted on the GWB at the low-frequency 
band, / > / crit ~ 2 x 10- 4 (O.lmm/I) 1 / 2 Hz [6] [7]. 

Despite the above interesting suggestions, however, quantitative prediction for 
the stochastic GWB in the brane-world cosmology is still under investigation. 
The spectrum of GWs generated by de Sitter inflation on a brane was first 
calculated by [4]. And later, the cosmological evolution of the GWs has been 
analytically studied in a very idealistic situation: the transition from the de- 
Sitter brane to the Minkowski brane [8] and the transition from the de-Sitter 
brane to the de-Sitter brane with different cosmological constant [9] . For more 
realistic situation of the cosmological evolution with equation of state p = wp, 
however, the wave equation of GWs cannot be separable and the analytical 
treatment is generally intractable. We must tackle the complicated partial 
differential equation directly. 

Theoretically, the evolution of GWs in the brane-world cosmology is expected 
to deviate from the standard four- dimensional theory in the following two 
aspects: 1) the non-standard cosmological expansion due to the bulk gravity 
and 2) the excitation of the Kaluza-Klein mode(KK-mode) of graviton. The 
former may enhance the amplitude of GWs and the latter may suppress or 
modulate the GW form on a brane. These effects are particularly significant in 
the high-energy regime of the universe. An important question is which effect 
is dominant and how the amplitude of GWs changes in a realistic situation of 
the cosmological brane-world. 

For this purpose, we solve the wave equation of GWs numerically in the 
Randall-Sundrum type single-brane model. It has been shown that the ex- 
citation of KK-modes are suppressed during inflation and zero-mode remains 
constant after inflation at super- horizon scales [4]. Thus we will consider the 
evolution of the GWs starting from an initially zero-mode state at super- 
horizon scales and focus on the behavior just after the horizon-crossing time. 
In the high-energy regime, the separation between zero-mode and KK-modes 
does not hold for the modes with wavelength shorter than the horizon. This 
implies that even in an initially zero-mode state, KK-modes are inevitably 
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excited when the perturbation crosses the horizon. Thus, we can observe how 
the excitation of the KK modes modifies the behavior of the perturbations. 

We set up the basic equations in section 2. Based on these equations, the nu- 
merical calculation of wave equation is performed and the results are presented 
in section 3. To understand the behavior on the brane, section 4 describes the 
semi-analytic treatment. Using the low-energy approximation, we derive an 
effective equation of the GWs on the brane. The resultant equation reduces 
to the ordinary differential equation and the numerical results of this equa- 
tion reasonably match those of the full numerical treatment. Final section 5 
is devoted to the summary and discussion. 



2 Basic equations 



In this letter, we specifically treat the single-brane model embedded in a five- 
dimensional anti-de Sitter space, in which the matter content on brane is 
simply given by a homogeneous and isotropic perfect fluid satisfying the equa- 
tion of state p = wp. Using the Gaussian normal coordinate to the brane, the 
background metric is given by 

ds 2 = -n 2 (t, y)dt 2 + a 2 (t, y)5 ij dx i dx j + dy 2 , (1) 

where the brane is located at y — 0. The lapse function n(t,y) and the scale 
factor a(t, y) for the background space-time are determined from the five- 
dimensional Einstein equation. In absence of the dark radiation, these quan- 
tities are written as follows [3] [4]: 

a(t,y) = ao(t) jcosh(^) - (l + sinh(^y) J , (2) 

a (t) (3) 

= e~ m + (2 + 3w)^Y smh(py), 
A 

where A > denotes the tension of the brane and p~ l represents the curvature 
scale of the anti-de Sitter bulk, which is related to the tension A by p = 
\J k\\/3. The scale factor ao(t) and the energy density on the brane p(t) are 
determined from the effective Friedmann equations [3]: 



H 2 = 



= Y p{t) f 1 + #) ' p = ~ 3H{1 + w)pit) ' (4) 
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The solution of the above equations is easily obtained and can be expressed 
in terms of the dimensionless variables r = t/t* and e = p/A: 



a (r) = a* (r 2 + 2cr - 2c) 3(1+w) , (5) 
p(r) =p* (t 2 + 2ct -2c)~\ (6) 

c^yi+f-i, (7) 

where the subscript * means the quantity evaluated at the time r = 1. Note 
that the Gaussian- normal coordinate has a coordinate singularity at y = y c , 
defined by a(r, y c (j)) = 0: 

^ c (r) = ilog{^(r + c) 2 }. (8) 



Now consider the tensor perturbation in the metric (1). The perturbed metric 
becomes 

ds 2 = -n 2 (t, y)dt 2 + a 2 (t, y)(5ij + E ij )dx i dx j + dy 2 . (9) 

The perturbed quantity satisfies the transverse and the traceless condi- 
tions, which is automatically gauge- invariant. Decomposing the perturbation 
in spatial Fourier modes as = E(t,y; k)e lkx iij, where is transverse- 
traceless polarization tensor, the wave equation for the Fourier component E 
becomes [4]: 

d 2 E fa n\dE n 2 fc2 ^ _ ^ f d 2 E fa' n'\dE) 
dt 2 V a nj dt a 2 " \ dy 2 \ a n J dy J 

where the dot and the prime respectively denote the derivative with respect 
to the time t and the bulk coordinate y. The above equation must be solved 
under the boundary condition at the brane, which follows from the perturbed 
Israel condition. Ignoring the anisotropic stress tensor ir^ arising from the 
matter on brane, we obtain the junction condition: 



dE 
dy 



= 0. (11) 

y=Q 



Hereafter, we will focus on the time evolution of the tensor perturbation in 
the radiation-dominated epoch(w = 1/3) and solve the wave equation (10) 
numerically. 
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3 Numerical analysis 



To solve the wave equation (10) numerically, we use the spectral method, 
which is a standard technique in the computational fluid dynamics [10] and 
recently becomes popular in subject of numerical relativity (e.g., [11]). To be 
precise, we adopt a Tchebychev collocation method with the Gauss-Lobatto 
collocation points, ye = cos(£n/N) with £ = 0, 1, • • • ,N. With this method, the 
quantity E(r,y) is first transformed to a set of variables E^{t) defined in the 
Tchebychev space, and the partial differential equation (10) can be regarded 
as the ordinary differential equations for Eg. Then, the time evolution of Eg, is 
followed by the Predictor- Corrector method based on the Adams-Bashforth- 
Moulton finite-difference scheme. 

As mentioned in section 1, the initial condition for the quantity E is set to 
the zero-mode solution in the inflationary era, i.e., E(ro,y) = constant. The 
calculations are then started before the wavelength of the zero-mode crosses 
the Hubble horizon. For convenience, we set the comoving wave number k to 
k = a^H*, that is, the G W mode just crosses the Hubble horizon when r = 1. 

The difficulty in present numerical calculation is that the computational do- 
main should be finite. In addition, the Gaussian normal coordinate (1) has a 
coordinate singularity (see Eq.(8)). Hence, we must introduce another bound- 
ary in the bulk and impose another boundary condition at this boundary 
for numerical purpose. Perhaps, the simplest choice is the static boundary, 
however, due to the narrow computational domain, the suppression of light 
KK-modes becomes significant and the artifitial reflection wave affects the 
physical brane. In this letter, to remedy this, we put a moving regulator brane 
at Vvcgij) = lUcij) inside the coordinate singularity and impose the Neumann 
boundary condition 1 : 

""a? < 12 » 

y=Vr C g (T) 

where n M is the normal vector to the boundary trajectory. To impose this 
boundary condition in the spectral method, we introduce the new coordinate 
z = 2y/y rcg (r) — 1, and solve (10) in the fixed range z : [—1, 1]. Thus, provided 
the location of the boundaries y = and y = y rC g(j), the only physical param- 
eter is e*, i.e., the dimensionless energy density p/X at the horizon-crossing 
time, which is directly related to the scale of the GW concerned. For large e*, 
the physical brane is very close to the coordinate singularity and the computa- 
tional domain of our code becomes narrow. This causes pathological behavior 



1 Restricting the parameter to e* < 1, the calculation with the static boundary can 
reproduce the results using the moving boundary at an early time, however, the 
late-time pahse of the calculation is significantly affected by the artifitial reflection 
wave. 
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near the coordinate singularity. Hence, in this letter, we restrict our analysis 
to the the small e* case. 
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Fig. 1. Numerical solution of GW on brane (solid), together with the 4-dimensional 
GW (long- dashed) and the result from the low-energy approximation (short- dashed, 
see Sec. 4). The upper right window is the zoom-up of the solutions in the range 
10 < r < 30. 

In Fig.l, the resultant waveform on the brane is depicted as function of time in 
the case with e* = 0.2, 7 = 0.8 and N = 64. Apart from the overall damping 
l/ao(t) due to the cosmo logical expansion, the amplitude of the numerical 
solution on the brane (solid) becomes smaller than that of the solution in 
the four- dimensional theory (long- dashed). This result simply reflects the fact 
that the localization of gravity is not fulfilled in presence of the p 2 -term (see 
Eq.(4)) and the GW can easily escape from the brane to the bulk. 

Fig. 2 shows the snapshots of the waveform plotted as function of the original 
coordinate y. In the low-energy limit (p/X <C 1) of the present case, the 
equation (10) becomes separable and the analytic solution can be obtained as 
E(r,y) = u(y)cf)(T). The mode function for the bulk coordinate u(y) is given 



by 



u (y) = constant 



(13) 



for the zero-mode and 




(14) 
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Fig. 2. Snapshot of the waveform as function of the bulk coordinate y in the case 
with the parameters (e*,7) = (0.2,0.8). 

for the KK-mode with the effective mass m. The function H^\z) denotes the 
Hankel function for the second kind. We attempted to approximate the wave 
form of the numerical result by the superposition of the mode funcitons (13) 
and (14). If we set the mass of the KK-mode to the Hubble parameter at the 
horizon crossing time, the mode functions roughly reproduce the waveform of 
the numerical result, although the agreement is at a qualitative level. It is clear 
from Figs. 1 and 2 that the created KK-modes is rather light, consistent with 
our numerical set-up with e* < 1, however, to obtain the KK-mode spectrum, 
a more quantitative analysis is needed. 

The introduction of the regulator brane in the bulk is certainly not a proper 
way to treat the boundary condition in the bulk. The regulator brane may 
eventually affect the physical brane located at y — 0. As shown in Fig.l, 
however, the excitation of the massive KK-modes becomes significant only 
at the horizon crossing time. The modification of the amplitude is almost 
determined around the horizon crossing and we need not follow the evolution 
of the perturbations for a long time. Thus, if the regulator brane is located 
sufficiently far away from the physical brane, we expect that we can safely 
neglect the effect of the regulator brane. To confirm this, in Fig. 3, we plot the 
damping factor, i.e., the ratio of the amplitude obtained from the numerical 
solution to that from the four- dimensional theory in the case of e* = 0.05, 
0.2 and 0.3. The number of collocation points used here is N = 64 except for 
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Fig. 3. Damping factor of the GW amplitude on brane as function of boundary 
parameter 7 

the cases of the boundary parameters 7 = 0.85 and 0.9 with iV = 128. Since 
the suppression of GW amplitude is effective only when the mode crosses the 
Hubble horizon, the damping factor evaluated at the end of the calculaiton 
is time-independent. The damping factor depends linearly on small e*, and 
tends to converge for 7 — ► 1. 



4 Result from low-energy expansion 

In order to understand the behavior of the numerical results in previous sec- 
tion, we employ a low energy expansion scheme to solve the wave equation 
approximately. In this treatment, the term p/X in effective Friedmann equa- 
tion (4) is assumed to be small. Thus, treating the dimensionless quantity 
e = p/X as small expansion parameter, we expand the tensor perturbation 
E(r,y), the scale factor ao( r ) and the Hubble parameter H(t): 



E(t, y) = E (t) + E 1 (r, y ) + E 2 (r,y) + --- , 



(15) 
(16) 
(17) 
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where E k (r,y), a^fY) and H k (r) are the functions of the order 0(e k ). Note 
further the fact that in the low energy approximation, the time-derivative 
of the lower-order quantities is comparable to the y-derivative of the higher- 
order ones [12]. This is because the time derivative is of the order of the Hubble 
parameter, while the y-derivative is of the order of the bulk scale /i. Thus, we 
have: 



-0(e). (18) 



Keeping the relation (18) in mind, one can expand the wave equation (10) 
in terms of e and find the higher-order solutions iteratively. Then, imposing 
boundary conditions, we derive an ordinary differential equations for E (t) 
[13]. 

The derivation of the equation for E (t) is straightforward but lengthly cal- 
culation. The final form of the equation up to the first order in e is 

□^0 = Ml g +W)2 (1 - 2B(y )) nlE 

- ^ {(2 + 3w) - (3w + 5)B(y )} n r E 

+ h§- Q -l(5 + 3w)(l + w)^(l-2B(y )) ] " * (19) 

k N 2 



/; 2 -(5 + 3w)(l + w)^(l- 2£?(y ))} H E 
+ {-2 ftl -3(l + ^(l-2%))UA| E , 




in the case of the static boundary. Here, we defined the operator O t and the 
quantity B as follows: 

S., (20) 



The variable yo denotes the position of the boundary determined at an initial 
time r = t , i.e., y = 72/c( r o)- If the regulator brane is moving, the equation 
contains the term proportional to y VC g(f) and it becomes rather complicated. 
Note also that this approximation will break down if the regulator brane be- 
comes far away from our brane: 

e ra rcg > p / A) (22) 

since the low energy approximation breaks down on the regulator brane. 

In spite of the above limitations of this approximation, the numerical results 
for Eq(t) quite well reproduces the full numerical simulation even in the case 
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of the moving boundary 2 (see short- dashed in Fig.l). Thus, at least in a 
qualitative level, we can use this approximation to understand the behavior 
of perturbations on the brane. 

The effective equation (19) contains a higher derivative term O 2 E which de- 
scribes the non-local effect caused by the propagation of the wave in the bulk. 
If we numerically solve the equation with this higher-derivative term, we find 
that the solution badly diverges in time. It might be caused by the trunca- 
tion of infinite numbers of higher derivative terms at a finite order. Hence 
we neglect the O 2 E term when we solve the ordinary differential equation 
(19) numerically. The first terms in the coefficients of HqEq and k 2 e~ 2ao Eo 
come from the modification of background Friedmann equation, that is, p 2 - 
term in Eq. (4). The second terms arise from the non-separable nature of the 
metric (1), which can be deduced from the fact that they vanish if we take 
w — — 1. These are the terms which excite KK-modes in the bulk and cause 
the dissipation of the perturbations on the brane. 

In the long wavelength limit k — > 0, the zero-mode solution E = const, 
is the solution for the effective equation (19). However, if the perturbation 
crosses the horizon, the zero-mode solution which satisfies O t Eq = cannot 
be a solution. Then the KK-modes are inevitably excited. The effects of the 
non-standard cosmological expansion tend to enhance the amplitude of the 
perturbation compared to the four-dimensional theory. On the other hand, 
both the numerical results and the low energy approximation show that the 
amplitude of the perturbations decreases, which implies that the influence of 
the excitation of KK-modes overcomes the effects of the non-standard cosmo- 
logical expansion. Therefore, the suppression of the amplitude of GWs can be 
understood as the consequence of the excitation of the KK-modes due to the 
non-separable nature of the bulk metric in the Gaussian normal coordinate 
defined with respect to the physical brane. 



5 Summary and Discussion 

In this letter, we have numerically investigated the evolution of GWs dur- 
ing the radiation dominated epoch in the Randall-Sundrum type single-brane 
model. Especially focusing on the behavior after the horizon-crossing time, we 
found that the amplitude of GWs on brane is suppressed compared to that 
of the four-dimensional theory. To interpret the numerical results, we also 
employ the low energy approximation and derive the effective wave equation 
on the brane (Eq.(19)). The solution of this equation reasonably agrees with 

2 Note also that the numerical integration of effective equation (19) reproduces the 
full numerical simulation with the static boundary if e* < 1. 
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numerical simulations and the suppression of GWs can be understood as a 
consequence of the excitation of KK-modes. Although the created KK-modes 
seem to be rather soft, the influence of the KK-modes on GWs overcomes the 
effect of the non-standard cosmological expansion arising from the p 2 -term. 
The suppression of GWs becomes significant as increasing the energy scales 
e*. Therefore, contrary to the four-dimensional prediction, the intensity of the 
stochastic GWB around the frequency / ~ / cr i t tends to decrease as increasing 
the frequency, which might provide an important clue to probe the presence 
or the absence of extra-dimensions. 

For more quantitative prediction for the spectrum of stochastic GWB, how- 
ever, the present numerical analysis has several limitations. While the evolved 
amplitude of GWs is examined at the relatively low energy scales e* < 1, 
corresponding to the low-frequency modes with / < f cr - lt , the high-frequency 
GWs with / > / crit are expected to suffer from the high-energy effects more 
seriously. As pointed out by [14], the created KK-modes does not simply es- 
cape away from the brane to the bulk in the high-energy region e* > 1. They 
could be turned back to the brane by the curvature scattering of the anti-de 
Sitter bulk, leading to the significant enhancement or the modulation of the 
GWs on brane. 

In order to investigate the degree of this effect precisely, the boundary con- 
dition imposed on the regulator brane might be inadequate. Instead of using 
the Neumann condition (12), a suitable choice of the boundary condition such 
as the non-reflecting boundary condition should be considered. Further, to 
impose a realistic boundary condition, the Poincare coordinate, which covers 
the wider region of the anti-de Sitter spacetime than the Gaussian normal co- 
ordinate, would be crucial in our numerical calculation. The implementation 
of these technical points is straightforward and the analysis is now in progress. 
We will report the results in a separate paper [13]. 
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